# -*- coding: utf-8 -*-
"""
Contains functions to test e-GARP and find the threshold e so that a dataset passes e-GARP.
The ``money'' version of functions takes as an input a N x N matrix called V where (i,k) 
reports the value of bundle C[i] at price level P[k]
"""


def garp( C,P,e = 1.0 ):
    '''Pass two N x L matrices. C is the consumption matrix and P is the price matrix.
    There are N observations and L goods. Will return True if the dataset satisfies
    e-GARP and False if it does not.'''
    V = C @ P.T
    return garp_money(V,e)
    
def garp_exact(C,P):
    '''Test for eact GARP. No CCEI.'''
    V = C @ P.T
    N = V.shape[0] #Number of obs
    
    #This loop will see if there are any cycles.
    for n in range(N):
        #Check if there is a cycle involving ob n.
        visited = N * [False] #If an element enters here we disregard it.
        queue = [n] #Observations left to check.
        while len(queue):
            m = queue.pop() #We will see what ob m is e-revealed preferred to
            for i in range(N):
                #The for loop does not consider obs we know aren't in a cycle.
                if (visited[i] == False) & ( i != m ): #Only check out obs that haven't been visited
                    if V[m,m] >= V[i,m]: #See if weak rp
                        #m is e-revealed preferred to i
                        if i == n: #If i is n then we've found a cycle unless the rp is weak.
                            if V[m,m] > V[n,m]:
                                return False
                        queue.append(i) #We know there's a path from n to i. 
                        visited[i] = True
    return True

def garp_exact_doc_cycle(C,P):
    '''Test for eact GARP. No CCEI. Returns a list of pairs of observations
    forming a cycle and the difference between the cost of a bundle and the dominated bundle.'''
    
    V = C @ P.T
    N = V.shape[0] #Number of obs
    
    #This loop will see if there are any cycles.
    for n in range(N):
        cycle_dic = { k: [] for k in range(N)} #Keep track of path from n to each ob k
        
        #Check if there is a cycle involving ob n.
        visited = N * [False] #If an element enters here we disregard it.
        queue = [n] #Observations left to check.
        while len(queue):
            m = queue.pop() #We will see what ob m is e-revealed preferred to
            for i in range(N):
                #The for loop does not consider obs we know aren't in a cycle.
                if (visited[i] == False) & ( i != m ): #Only check out obs that haven't been visited
                    if V[m,m] >= V[i,m]: #See if weak rp
                        #Document this
                        cycle_dic[ i ] = cycle_dic[m] + [[ m, i, V[m,m] - V[i,m] ]]    
                        #m is e-revealed preferred to i
                        if i == n: #If i is n then we've found a cycle unless the rp is weak.
                            if V[m,m] > V[n,m]:
                                return cycle_dic[n]
                        queue.append(i) #We know there's a path from n to i. 
                        visited[i] = True
    return {}

def garp_money( V, e = 1.0):
    '''
    Pass an N x N matrix called V where (i,k) reports the value of bundle C[i] 
    at price level P[k]. 
    
    e is a number between 0 and 1. It is the efficiency level. 
    Note that if you pass e==1 it may still be possible that you don't pass GARP.'
    
    Returns True if the dataset passes GARP at efficiency level e. Otherwise False.
    '''
    
    N = V.shape[0] #Number of obs
        
    #This loop will see if there are any cycles.
    for n in range(N):
        #Check if there is a cycle involving ob n.
        visited = N * [False] #If an element enters here we disregard it.
        queue = [n] #Observations left to check.
        while len(queue):
            m = queue.pop() #We will see what ob m is e-revealed preferred to
            for i in range(n,N):
                #The for loop does not consider obs we know aren't in a cycle.
                if visited[i] == False: #Only check out obs that haven't been visited
                    if e * V[m,m] > V[i,m]: #See if rp
                        #m is e-revealed preferred to i
                        if i == n: #If i is n then we've found a cycle
                            return False
                        queue.append(i) #We know there's a path from obnum to i. 
                        visited[i] = True
    return True

def garp_find_e( C,P ):
    '''C and P are N x L matrices where N is the number of obs and L is the 
    number of goods. C represents the consumption bundles and P the prices.
    Finds the efficiency level at which the data just passes e-GARP.'''
    V = C @ P.T
    return garp_money_find_e(V)
    
    
def garp_money_find_e( V ):
    '''V is an N x N matrices where N is the number of obs. V[n,m] represents
    the value of bundle n evaluated at price level m. 
    Finds the efficiency level at which the data just passes e-GARP'''
    
    high_e = 1.0
    low_e = 0.0
    e = 1.0
    #Will binary search until we have narrowed the window enough.
    while high_e - low_e >= 2**(-15):
        if garp_money( V, e ):
            low_e = e #Our e passes so raise the lower bound.
        else:
            high_e = e #Our e failed so lower the upper bound.
        e = (high_e + low_e) / 2.0 #New e is the midpoint.
    return e




